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ABSTRACT 


A new method to compute lunisolar perturbations in satellite motion is proposed. 
The disturbing function is expressed by the orbital elements of the satellite and the 
geocentric polar coordinates of the moon and the sun. The secular and long-periodic 
perturbations are derived by numerical integrations, and the short-periodic perturba- 
tions are derived analytically. The perturbations due to the tides' can be included in 
the same way. 

In the Appendix, the motion of the orbital plane for a synchronous satellite is 
discussed; it is concluded that the inclination cannot stay below 7°. 


RESUME 

On propose une nouvelle mgthode de calcul des perturbations luni-solaires dans 
le mouvement d'un satellite. La fonction perturbante est exprimge en fonction des 
elements orbitaux du satellite et des coordonnges polaires ggocentriques de la 
lune et du soleil. On deduit les perturbations s^culaires et a longue pdriode par 
integrations numdriques et les perturbations a courte p§riode par 1'analyse. Les 
perturbations dQes aux mardes peuvent £tre incluses de la m§me manifere. 

On considfcre, dans 1'Appendice, le mouvement du plan orbital pour un satellite 
synchrone; on conclut que 1 1 incl inaison ne peut rester infgrieure a 7°. 
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KOHCFlEKT 

npeflnaraeTCH HOBbiw MeTOfl ajih BMHMCJieHWH jiyHHO-cojiHenHbix 
B03MymeHMH. 3>yHKu;MH bo 3 MymeHMH BbipaJKaeTCH nyTeM opdHTajibHbix ajieMeH- 
tob cnyTHMKa m reopeHTpwiiecKMX KOopflMHaT jiyHbi n cojiHua. BeKOBbie m 
flOJiro-nepwoflHqecKMe bo 3 Mym eHnn bmbohhtch c noMombio HMCJieHHbix mht e— 
rpMpoBaHHM, a KopoTKO-nepMonymecKMe bmbohhtch aHaxiwTpmecKii . Bo3iwy- 
UieHHH BbI3bIBatOTCH npHJIWBaiWH,-'-M '.MO ryT .-6bl Tb . BKJIfOUeHbl' . no B.O 6HbIM o6pa30M. 

B npMjioaceHHM oScyxaaeTCH HBHxeHHe opCnTanbHOM naocKOCTn am 
CMHxpoHHoro cnyTHMKa; npnxonnTca k BbiBoay hto h£k/ioh He MOsceT 

n O 

ocTasaTbCH HMSce / 
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A NEW METHOD TO COMPUTE LUNISOLAR PERTURBATIONS 
IN SATELLITE MOTIONS 

Yoshihide Kozai 
1. INTRODUCTION : 

In an earlier paper (Kozai, 1959), the principal secular and long-periodic terms 
of the disturbing function for the satellite motion due to the lunisolar gravitational 
attractions were given as a function of the orbital elements of the satellite, the sun, 
and the moon. The parallactic term in the disturbing function was later included by 
Musen, Bailie, and Upton (1961) . 

The disturbing factors are 1. 63 X 10 5 /n 2 and 0. 75 X 10 5 /n 2 for the moon and 
the sun, respectively; when the mean motion n of the satellite is expressed in revolu- 
tions per day, some of the amplitudes for the perturbation terms may become rather 
large, and they must be computed to 3 or 4 figures even if the satellite is close to the 
earth. 

However, in order not to make the expressions too complicated, it is usually 
assumed that the orbital elements of the sun and the moon in the disturbing function 
are constant except for the lunar ascending node and perigee, which are assumed to 
move as linear functions of time . 

It is a good assumption for computing the solar perturbations. However, that is 
not the case for the lunar perturbations, because it is known that the lunar eccentricity 
and the lunar inclination to the ecliptic change between 0. 045 and 0. 065 and between 
5° O'. 5 and 5° 17', respectively, and there are some large inequalities with amplitudes 
as large as a few tenths of a degree in the lunar mean longitude. 


This work was supported in part by grant NGR 09-015-002 from the National 
Aeronautics and Space Administration. 
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Therefore, many celestial mechanicians conclude that they must depend on numer- 
ical integration methods if they need to compute precisely the lunisolar perturbations, 
particularly when the satellite is high. 

In this paper, I propose a new method to compute the lunisolar perturbations. The 
disturbing function is expressed as a function of the orbital elements of the satellite 
and the polar coordinates of the sun and the moon; the latter can be found in any 
astronomical almanac. From the disturbing function of this form, the short-periodic 
terms can be eliminated analytically by taking their average with respect to the mean 
anomaly of the satellite. 

Then the equations of variations of the orbital elements must be integrated 
numerically to derive long-periodic perturbations. However, since the short-periodic 
terms have been eliminated, a half-day or even a day can be adopted as a step of the 
integration. 

This method can be applied to compute the lunisolar perturbations for high satel- 
lites such as a synchronous one, and the short-periodic perturbations can be derived 
analytically. 

The perturbations due to the tidal deformations of the earth can be included in the 
same way; their expressions are also given. 
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2. DISTURBING FUNCTION 


The expression of the disturbing function due to the lunisolar gravitational attrac- 
tions for the satellite motion can be derived in exactly the same way as that used to 
derive the disturbing function for the lunar motion due to the sun. 


The equations of motion for the satellite under the attractions of the spherical 
earth and the spherical moon are written as 


d 2 r 
dt 2 


= G — 
9r 


2 3 

™ p 2 ( cos s ) + m ' P 3 ( cos s ) + ' • ’ 

r' r' 


( 1 ) 


where r and r are the geocentric radius vector and distance of the satellite, r' is the 
geocentric distance of the moon, m and m' are the masses of the earth and the moon, 
G is the gravitational constant, cos S is the direction cosine between the radius to the 
satellite and that to the moon, and P 2 , P^, • • • are Legendre polynomials. 


The first term in parentheses in (1) expresses the potential of the spherical earth; 
when the other terms are neglected, the equation becomes that for the two-body 
problem. The second term in parentheses gives the principal part of the disturbing 
function. 


The disturbing factor, which is the ratio of the second term to the first, is 
expressed as 



where a and a' are the semimajor axes of the satellite and the moon, respectively. 
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By use of Kepler's third law, 


2 3 „ 

n a = Gm , 


n ,2 a' 3 = G(m + m') , 


( 3 ) 


expression (2) becomes 


m + 


1 — (— ) 2 
m' Vn / 


( 4 ) 


The ratio of the third term to the second is roughly equal to a/a' and is not a very 
small quantity. 

The disturbing function due to the sun and the disturbing factor for the sun are 
derived in the same way. However, as Kepler's third law for the sun is written as 


n' 2 a r 3 = Gm' , 


( 5 ) 


the disturbing factor becomes 


(n'/nr 


( 6 ) 


And the ratio of the third term to the second is extremely small. 

2 3 

When Gm is replaced by n a , the disturbing function can be written as 


R = 


, 2 _ 2 / a '' 

» P r [p- 


> \ 3 r 


P 2 (cos S) + p- P 3 (cos S) + 


/ \2 

( r t 

\r' ) 


P 4 (cos S) + 


( 7 ) 


where 


P = 


m 


m + m' 


for the moon , 


P = 1 


for the sun 


( 8 . 


4 



The geocentric rectangular coordinates of the moon or the sun are expressed with 
the geocentric distance r' , right ascension a, and declination 6 as 

x' = r' cos 6 cos a , 
y' = r' cos 6 sin a , 

z' = r' sin 6 , (9) 

whereas the coordinates of the satellite are expressed with the geocentric distance r, 
the argument of latitude L, the longitude of the ascending node £2, and the orbital 
inclination i as 


x = r cos L cos £2 - r cos i sin L sin £2 , 
y = r cos L sin £2 + r cos i sin L cos £2 , 

z = r sin L sin i (10) 

Then the direction cosine is expressed as 

cos S = (xx' + yy' + zz')/rr' 

= cos 6 cos (£2 - a) cos L 
+ [- cos 6 cos i sin (£2 - a) + sin 6 sin i] sin L 
= A cos L + B sin L , (11) 

where 

A = cos 6 cos (£2 - a) , 

B = - cos 6 cos i sin (£2 - a) + sin 6 sin i . (12) 
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By use of the last expression for cos S in (11), the following formulas are derived: 


9 a 2 + r^ a 2 - R^ 

cos^ S = + g cos 2L + AB sin 2L , 

cos 2 S = | (A 2 + B 2 )(A cos L + B sin L) 

A 2 - 3B 2 3A 2 - B 2 

+ A cos 3L + B sm 3L , 

4 4 ’ 

cos 4 S = | (A 2 + B 2 ) 2 + (A 2 + B 2 ) [(A 2 - B 2 ) cos 2L + 2AB sin 2L] 

O Z 

. A 4 - 6A 2 B 2 + B 4 A 2 - B 2 . . 

+ cos 4L + - AB sm 4 L , 

O 6 

cos 5 S = 4 (A 2 + B 2 )(A cos L + B sin L) 

O 


5 sin sin 

In the expression of cos S, terms of 3L and 5L are neglected. 
^ ’ cos cos & 


Then the expressions for Legendre polynomials are derived as 

P 2 (cos S) = | [3(A 2 + B 2 ) - 2] +| [(A 2 - B 2 ) cos 2L + 2AB sin 2L] , 

P Q (cos S) = § [5 (A 2 + B 2 ) - 4] (A cos L + B sin L) 

o o 

+ | [(A 2 - 3B 2 ) A cos 3L + (3A 2 - B 2 ) B sin 3L] , 

O 

P 4 (cos S) = ££ [35(A 2 + B 2 ) 2 - 40(A 2 + B 2 ) + 8] 

+ ~ [7(A 2 + B 2 ) - 6] [(A 2 - B 2 ) cos 2L + 2AB sin 2L] 
lb 

+ || [(A 4 - 6A 2 B 2 + B 4 ) cos 4L + 4(A 2 - B 2 ) AB sin 4L] , 

Pj_ (cos S) = || [21(A 2 + B 2 ) 2 - 28(A 2 + B 2 ) + 8] (A cos L + B sin L) . (14) 

Here also, terms of 3L and 5L are neglected in the expression of Pg(cos S). 


Now the disturbing function is expressed by the orbital elements of the satellite 
and the polar coordinates of the moon and the sun. 
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3. SECULAR AND LONG-PERIODIC TERMS 


Secular and long-periodic terms of the disturbing function are derived by taking 
their average with respect to the mean anomaly A by the following formula: 


2ir 

h J (!) cos 1,11 ■ il 

0 


, ,m (n + 2)(n + 3)(n + m + 1) 
ml 


F 


/ m - n 

V 2 



m + 1 



(15) 


where f and e are the true anomaly and the eccentricity, respectively, and F is the 
hypergeometric function. 


2 4 3 3 3 5 4 

If (a/a' ) e , (a/a' ) e , (a/a' ) e , and (a/a' ) are neglected, the secular and 

long-periodic terms of the disturbing function are written as 



= | [3(A 2 + B 2 ) -2] (2 + 3e 2 ) 
- B 2 ) cos 2 oj + 2AB sin 2u] 


~(^T e (r r )(a 7 ")^ 3 ^ 4 + 3e2 H 5 ( A2 + b2 ) “ ( A cos w + B sin w) 

+ 35e 2 [(A 2 - 3B 2 ) A cos 3u + (3A 2 - B 2 ) B sin 3co]} 

+ ^(fr) 2 (^) 2 {[35(A 2 + B 2 ) 2 - 40 (A 2 + B 2 ) + 8] (1 + 5 e 2 ) 

+ 35e 2 [7(A 2 + B 2 ) - 6 ] [(A 2 - B 2 ) cos 2co + 2AB sin 2u>]} 

- HI e (~ 7 ~j (jf) [21(A 2 + B 2 ) 2 - 28(A 2 + B 2 ) + 8] (A cos w + B sin w ) , 

(16) 


where w is the argument of perigee for the satellite. 


7 



By use of Lagrange's equations, 


da 

_ _2_ 

3R 




dt 

na 

3 St 

y 



de 

_ 

2 

e 

3R 


2 

- e 

dt 


2 

d£ 


2 


na 

e 


na 

. e 

di 


cos i 


3R 

dt 

2 


2 


3co 


na 


dco 

dt 


cos 1 


3H , Vl - e 2 3R 


2 2 . . 
na VI - e sin i 


3i 


2 3e ’ 
na e 


(in 

dt 


2 \f. 2 . . 

na VI - e sm i 


3R 

3i ’ 


1 - e* 3R 2 3R 


di = 

dt n 2 3e na 3a ’ 
na e 


the following equations are derived: 


^ = 0 
dt ’ 


3R 

dn ’ 


de 

dt 


n' 2. , , ,,3 _ J ~ I 2 

a (a' /r ) 6 V 1 - e 

n J 

^ e[(A 2 - B 2 ) sin 2<o - 2AB cos 2to] 


- ^(r 7 ") ("a 7 ") + 3e2 )^ 5 ( A2 + B 2 ) “ 4] (A sin to - B cos to) 

+ 105e 2 [(A 2 - 3B 2 ) A sin 3<o - (3A 2 - B 2 ) B cos 3to]} 


(17) 


(18) 


8 



f \ 2 


W (p) (f) e[7 ( A2 + ® 2 ) ■ 6] [ ( a2 - b2 ) sin 2w - 2AB cos 2«1 

3 3 

iH (f) (f) [21(a2 + ® 2)2 " 28(a2 + b2) + 81 (A sin <*> - B cos <*>) > 

(19) 


di , e A . de 

37 + o cot i -t— 

dt . 2 dt 

1 - e 


r ,2 . 3 

- V(f) p 


yl-e sini 


! ( A i +B i^ < 2 + 3e2 > 


15 2 

T e 


. 9A _ 9B \ . ^ / . 9B _ 9A , . „ 

A FT B ar 7 cos2w + A an + B a^ |sm2w 


, 5 /a' 

64 6 Vr' 


3(4 + 3e ) 


(15A 2 + 5B 2 - 4)f| + 10 Ab|| 


COS co 


(5A 2 + 15B 2 - 4) f + 10AB f| 


sin co 


+ 105e 


(A 2 r 2 \ ^ A - 2AR ^ B 

(A B > an 2AB an 


cos 3 co 


( A 2 - B 2 ) || + 2AB || 


sin 3 co 


- S (f ) 2 (f ) 2 ( 2[7 < a2 + ® 2) - 41 ( A 


b|§ ) (1 + 50 


+ 17e <2 


_ , .3 9A 3 9B 

7,A an“ B 0^ 


o [ a — R — 

3 1 A an B an 


cos 2 


CJ 


7 | qa 2 R M. + A 3 + q AR 2 M + rS 9A 

7 l3A B an + A an + 3AB an + B an 


6 (A i +B i 


sin 2co; 
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+ ii e iff (r) 3 I 28[3 < a2 + b2 > - 21 ( A m + B H) < A c ° s w + b sin u) 


+ [21(A 2 + B 2 ) 2 - 28(A 2 + B 2 ) + 8] ( |~ cos w+ sin to) 


(20) 


d£2 
dt I 


n' 2 / a' \ 3 1 1 9B 

n \ r' / " / r sin i 3i 


VTT7 


= | B(2 + 3e 2 ) - ^ e 2 (B cos 2 to - A sin 2oo) 


|| e fp-) (jfj {(4 + 3e 2 )[10AB cos to + (5A 2 + 15B 2 - 4) sin to) 


- 35e iS [2AB cos 3 go - (A 2 - B 2 ) sin 3 to]} 



- 7e 2 [2B(7B 2 - 3) cos 2co - (7A 2 + 21B 2 - 6) A sin 2 co] } 
ill 6 (p) 3 (^-) 3 {28AB[3(A 2 +B 2 ) -2] cos 2 to 


+ (21A* + 126A 2 B 2 + 1056* - 28A^ - 84B^ + 8) sin 2co} , (21) 


du , ■ . dfi' 

dt" + 008 dt" 


/«'\ 3 


P Vi- 


n \r 


| [3(A 2 + B 2 ) - 2] + -p [(A 2 - B 2 ) cos 2co + 2AB sin 2to] 


|f (p~) (fr) | {(4 + 9e 2 )[5(A 2 + B 2 ) - 4] (A cos to + B sin to) 


+ 35e 2 [(A 2 - 3B 2 ) A cos 3u + (3A 2 - B 2 ) B sin 3 to]} 
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where 


+ 7 [7 (A 2 + B 2 ) - 6] [(A a - B^) cos 2u + 2AB sin 2 w j} 


l2 „2. 



( 22 ) 


di 

dt 


" ,2 , , » 3 

V F f 


- \ (7 + 3e 2 )[3(A 2 + B 2 ) - 2] 


- (1 + e 2 )[(A 2 - B 2 ) cos 2 co + 2AB sin 2co] 


+ If f-T-) ~ {(4 + 29e 2 )[5(A 2 + B 2 ) - 4] (A cos co + B sin co) 
v r / v a / g 


+ 35e 2 [(A 2 - 3B 2 ) A cos 3co + (3A 2 - B 2 ) B sin 3co]} 


- 4 (fr) (f-) (3(3 + 5e 2 )(35(A 2 + B 2 ) 2 - 40(A 2 + B 2 ) + 8] 

+ 35e 2 [7 (A 2 + B 2 ) - 6] [(A 2 - B 2 ) cos 2co + 2AB sin 2co]} 

3 3 

+ •*§§(£■) ^ [21(A 2 + B 2 ) 2 - 28(A 2 + B 2 ) + 8] (A cos co + B sin co) , 


(23) 


3A 

dn 

8B 

dn 

9R 

at 


- cos 6 sin (£2 - a) , 

- cos 6 cos i cos (£2 - a) , 

cos 6 sin i sin (£2 - a) + sin 6 cos i . 


(24) 
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These equations must be integrated numerically together with other terms, par- 
ticularly the secular terms due to J , because the secondary effects in the secular 
motions due to the changes of the orbital elements caused by the lunisolar perturbations 
may become as large as the direct effects. 

For close satellites, it is not necessary to include all the terms in equations 

2 3 

(19) through (23). In fact, since Gm has been replaced by n a , the expressions on 

2 

the right-hand side of the equations are accurate up to J^/ a . 

However, for high satellites such as a synchronous one, for which a/ a' is 0. 11 for 
2 -5 

the moon, and J /a is 2. 5 X 10 , all the terms must be included. For this case, a 

A 

more accurate expression for Kepler's third law for the moon, 

1.002723 n /2 a' 3 = G(m + m # ) , (25) 

must be used. If the value 1/81. 302 is adopted for m' /m, p for the moon takes the 
following value: 


p = 0.0121835 . 


(26) 



4. EXPRESSIONS FOR SMALL ECCENTRICITY AND INCLINATION 


As most of the existant and planned synchronous satellites have orbits that are 
nearly circular and nearly equatorial, expressions (19) through (23) are simplified by 
assuming that e and i are very small, although the inclination cannot remain very small, 
as explained in the Appendix. 

When sin 2 i, e 2 sin 2 i, (a/a') e sin 2 i, (a/a') e 2 , (a/a ') 2 sin 2 i, (a/a ') 2 e 2 , 

3 3 3 

(a/ a') e , and (a/a' ) e sin i are neglected, expression (16) becomes 


„ ,2 2 

R = n' a 




[3 cos 2 6 (2 - sin 2 i) + 6 sin 2 i sin 2 6-4 


2 2 2 
- 6 sin i sin 26 sin (n - a) + 3 sin i cos 6 cos 2 (tt - a)] • (2 + 3e ) 


15 2 2 

+ -jr e [cos 6 cos 2 (to + ft - a) + sin i sin 26 sin (2to + ft - a)] 

O 


15 

32 


/ sf \ / 3 . \ 2 

( —) ( — e [2 cos 6 (5 cos 6-4) cos (to + ft - a) 
y r / \ ci j 


+ 4 sin i sin 6 (5 cos 6-2) sin to- 5 sin i sin 26 cos 6 sin (to +2 ft - 2a)] 


2 2 

+ — (-^ 7 ) \ t) [35 cos^ 6-40 cos 2 6 + 8-10 sin i sin 26 (7 cos 2 6-4) 
V3. j \ r / 


64 \z 


sin (f2 - a)] 


105 

128 


/a\ 2 /a'\ 2 4 2 I 

( — J ( — ) (21 cos 6-28 cos 6 + 8 ) cos 5 cos (to + f2 - a) 

\ 3 . / y r / 


(27) 
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Then the equations of variations of the orbital elements are derived with sufficient 
accuracy as follows: 


de n' 2 / a' \ 3 2 f 15 , 2 _ . _ 

dt = ~ W) f 3 * 1_e e[cos 5 sin 2 


- sin i sin 26 cos (2J2 + f2 - a)] 


(to + £2 - a) 


- ~ \~t) f^r) [2 cos 6 (5 cos 2 6-4) sin (co + J2 - a) 
oZ \3. / \JC J 


4 sin i sin 6 (5 cos 2 6-2) cos to + 5 sin i sin 26 cos 6 cos (go + 2f2 - 2a)] 


105 

128 


3 /„/\ 3 


/ a \ /a' \ 

\a'/ \r' / 


(21 cos 6-28 cos 6 + 8) cos 6 sin (go + £2 - a)j> , 

(28) 


di 

dt 


2 3 ( 

| (~r) P I ^ [sin 26 cos (S2 - a) + sin i cos 2 6 sin 2 (£2 - a)] 

Vl - e 2 * 


15 / a \ 2 /a' \ 2 2 

+ ( — ) ( — ) sin 26 (7 cos 6-4) cos (t2 - a) 

32 \ a / \ r / 


dt2 

dt 




2 2 
T [- cos i cos 6 + 2 cos i sin 6 


(29) 


- cot i sin 26 sin (£2 - a) + cos i cos 6 cos 2(0 — a)] 

- -|| (y-) cot i ' sin 26 (7 cos 2 6-4) sin (f2 - a) ( , 


(30) 


lOJ 


dS2 


cos 1 


dt dt 

2 /_/\ 3 


+ P V 1 " ^ | f t 3 cos 2 6-2-3 sin i sin 26 sin (£2 - a)] 

15 2 

+ — [cos 6 cos 2 (go + S2 - a) + sin i sin 26 sin (2co + £2 - a)] 
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di 

dt 


15 

32 


(^ 7 ) (-^t) — [2 cos 6 (5 COS 2 6-4) cos (to + O - a) 
\ r / V 3. / g 


+ 4 sin i sin 6 (5 cos 6-2) sin to - 5 sin i sin 26 cos 6 sin (to + 20 - 2a)] 
- ■jH' (yaJ~) e ^ 21 cos4 & - 28 cos 2 6+8) cos 6 cos (to + O - a)] j , 


(31) 


n' 2 / a' \ 2 f 7 2 

— ( p-l (3 < - — [3 cos 6-2-3 sin i sin 26 sin (O - a)] 

15 2 

- — [cos 6 cos 2 (to + O - a) + sin i sin 26 sin (2 to + O - a)] 

15 1 / a \ ( a \ r 0 c /f - 2 c .. . . 

+ — — ( — ) ( — J 12 cos 6 (5 cos 6-4) cos (cj + ft - a) 

g \a / \ r / 

+ 4 sin i sin 6 (5 cos 2 6-2) sin to - 5 sin i sin 26 cos 6 sin (to + 20 - 2a)] 

o /a\ 2 /a' \ 2 4 2 

- 1 ( -4-J ( =r) (35 cos 6-40 cos 6 + 8) 

8 \a / \r / 7 


105 

128 


3 3 \ 

( -V) — (21 cos 4 6-28 cos 2 6 + 8) cos 6 cos (to + O - a)J- . 

\ a / \r / e l 


(32) 
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5. SHORT-PERIODIC PERTURBATIONS 


For some satellites, especially high ones, it is necessary to include, for precise 
orbit computations, short -periodic perturbations due to the lunisolar attractions. 


If only Pg (cos S) is included, the short -periodic term of the disturbing function 
is written as 


V "' 2 * 2 (£)% | f [3<A 2 + B 2 ) - 2] [(f) 2 -l-| 


3 2 

~ e 


+ | (A 2 - B 2 ) 


cos 2L - j e 2 cos 2 


+ f AB 


^ ) sin 2L - ^ e 2 sin 2 


(33) 


When the relations 


/r \ 2 3 2 e 2 

J = 1 + g e - 2 e cos 2 - — cos 2 2 + •■ ■ , 


( r \ 2 5 2 / 5 o\ 2 

— J cos 2f = g e - 3e cos i + ( 1 - — e 1 cos 2 f + e cos 3i + e cos 4i + ■ • • , 

2 

sin 2f = - 3e sin 2 + ^1 - e 2 ^ sin 2 2 + e sin 32 + e 2 sin 4 2 + •• • , 


(34) 


2 

are used and when e is neglected, the expression of R^ becomes 
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, 2 2 / a' \ 

V a (—) 


i 2 + B 2 ) cos St 


+ ^ (A 2 - B 2 ) [cos 2 (St + oj) - 3e cos (St + 2co) + e cos (3St + 2oo)] 

3 ) 

+ - AB [sin 2 (i + co) - 3e sin ( St + 2co) + e sin (3i + 2u>)] )■ . (35) 

2 j 

Then Lagrange's equations are easily derived. When the right-hand sides of the 
equations are integrated, it is assumed that all the variables except the mean anomaly 
St of the satellite are constant. By this approximation, n'/n is neglected in the results. 
Also, terms with e as factors are neglected. 


The results are as follows: 


9 Q 

rlo Q / n' \ / \ O O 

a = 2 ( n ) (f) P [(A " B ) COS 2 + “> + 2AB sin 2 (Jf + «)] , (36) 


de = ( n“) 2 (F) 3p | 1-|( a2 + B 2 ) cosi 
+ - (A 2 - B 2 )[9 cos (i + 2w) + cos (3i + 2u)] 


+ ^ AB [9 sin (i + 2w) + sin (3St + 2co)]> , 


edi = (f ) 2 (f ) 3(3 I" 1-|( a2 + B 2 ) Sini 
+ ^ (A 2 - B 2 ) [9 sin (l + 2co) - sin (3f + 2co)] 

- ^ AB[9 cos (St + 2cu) - cos (3f + 2u)l | , 
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2 3 ( 

d(i + q) = (p-) (p) p \-f [(A 2 - B 2 ) sin 2 (i + to) - 2AB cos 2 (i + to)] 

2 f 

+ — cot i gj- [B sin 2 (i + to) + A cos 2 (i + to)] > , (39) 

2 3 

dn = - | (p) (p) P -pr If [B sin 2 (i + to) + A cos 2 (i + to)] , (40) 

di = | (p) (p~) P ppy | cos i [(A 2 - B 2 ) cos 2 (i + to) + 2AB sin 2 (i + to)] 

- ( A l§- B li) sin2 < < + “> + ( A fi + B M) cos 2 (f + “'} ■ 

(41) 


To derive expression (39) for d(i + to), the fact must be taken into account that the 

3 del 

mean motion changes by - — n — , which is computed from (36) . 

u cL 


Expressions (36) through (39) can be combined to express the perturbations in the 
radius dr, in the argument of latitude dL, in the radial component of the velocity dr, 
and in the time derivative of the true anomaly df, by the following relations: 


dr = — da - a cos f de + 
a 


ae sin f ^ 

VT7 , 


y 


dL = dco + ( I* + — ^~2 
1 - e 


sin f • de + 




di 


df = - | ^ Vl - e 2 da - 2ne sin f (p) di 


+ n (—) (2 cos f + e cos 2f) de , 

VT7 
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sin f da + ane 


cos f • di 



( 42 ) 


If the eccentricity is neglected, the results are expressed as follows: 

^ = - (^-) 2 (|^-) 3 (3 1 - | (A 2 + B 2 ) + (A 2 - B 2 ) cos 2 (i + w ) 

+ 2AB sin 2 (I + go) , (43) 

2 3 

dL = -i (p-) (3 | 11(A 2 - B 2 ) sin 2 (i + co) - 22AB cos 2 (i + co) 

+ 6 cot i ~ [B sin 2 (f + co) + A cos 2 (i + co] j , (44) 

2 3 

dr = 2na (^-) (^-) (3 [(A 2 - B 2 ) sin 2 (i + u>) - 2AB cos 2 (i + co)] , 

(45) 

df = 4 n (if) (F-) 3 f 3 ^ 2 - 3 ( a2 + b2 )J 

+ 11 (A 2 - B 2 ) cos 2 (i + co) + 22AB sin 2 (i + co)} . (46) 

If the inclination is also small, these results can be written as 
2 3 

P 1 - 7 > cos 2 6 + cos 2 5 cos 2 (L + £7 - a)l , (47) 
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( 48 ) 


2 3 

dL = | (p) (p-) (3 [11 cos 2 6 sin 2 (L + f2 - a) 


+ 3 cot i sin 26 cos 2 (L + n - a)] , 

,2 ,3 

dr = 2na {~) (p-) (3 cos^ 6 sin 2 (L + 12 - a) , 

,23 

df = | n (p) (p-) p [4(2 - 3 cos 2 6) + 11 cos 2 6 cos 2 (L + n - a)] 


Also for the small-inclination case, expressions (40) and (41) are written as 

dn = -f © 2 (It ) 3 f 5 cos 2 < L + n - “> > 


di 


i = |- (p - ) P sin 25 sin (2L + n - a) 


,2 ,3 

' \ , n ' « 


(49) 

(50) 


(51) 

(52) 
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6. PERTURBATIONS DUE TO THE TIDAL DEFORMATIONS 


Besides the direct lunisolar perturbations, the motion of the satellite is disturbed 
by the tidal deformations of the earth due to the attractions of the moon and the sun 
(Kozai, 1965). 

Owing to the tidal potential 

2 

a 

Gm' — P cos S) , (53) 

r' 

where a g is the earth's radius, the geoid is deformed by 
4 

k 2 m" "7% P 2 (COS S) ’ (54) 

r' 

where is Love's number. And owing to the friction, the rotation carries the tidal 
bulge forward; the tide is high not when the moon is directly overhead but at some time 
delay At. 


The geopotential is also changed, by 


Gm/ 


S-T k 2 P 2 (C ° S S) 


(55) 


Therefore, the disturbing function due to the tidal deformation is written as 

3 

R = n' 2 pr 2 (fr) k 2 P 2 (cos S) . (56) 

However, owing to the time lag of the tide, the right ascensions of the sun and the 
moon should be shifted by n_ • At, where n ffi is the rotational angular velocity of the 

T37 w 

earth. 
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When this effect is taken into account, cos S can be written as 

* * 

A = cos 6 cos (S2 - a ) , 

* * 

B = - cos 6 cos i sin (J2 - a ) + sin i • sin 6 , (57) 


where 


a = a + n^ • At . 


(58) 


As the relations 


2tt 


eJ (?) 


di = (1 - e ) 


V "3/2 


2tt 


eJ (?)' 

0 


cos 2f ■ di = 0 


(59) 


hold, then after the short-periodic terms are eliminated, the disturbing function is 
written as 

5 


(f) 3 W-*v 3/2 

a 


| (A* 2 + B* 2 ) - \ 


(60) 


Lagrange's equations for this case are then easily derived as 


^=0 
dt ’ 


dt u ’ 


di _ 3 _J_ ( * 9A* * 9B* 

dt ~ 2 sin i \ d£2 9f2 , 
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5 


( 61 ) 


( 62 ) 
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APPENDIX 


MOTION FOR THE ORBITAL PLANE OF A SYNCHRONOUS SATELLITE 


As a synchronous satellite has the large value of 6. 65 equatorial radii for the 
semimajor axis, the disturbing factors due to 3^, the moon, and the sun are of the 
same order of magnitude: 

J 2 /a 2 = 2.46X 1(T 5 , 

(3(n' /n) 2 = 1. 63 X 10 5 for the moon , 

(3(n' /n) 2 = 0. 75 X 10 3 for the sun . (63) 

When both the eccentricity and the inclination are very small for the satellite, 
the main term of the disturbing function is written as (Kozai, 1959) 


„ 22 2 / 1 3 . 2 A / 3 2 

R = n a \ 2~2 Sm V [ 1 + 2 6 


n' 2 a 2 p(l4e 2 | (1-e' 2 )- 3 / 2 


1 / ■ 3 . 2 A / 3.2., 

4 l 1 “2 Sm 1 l 1 “ 2 Sin 1 


3 3 2 2 

+ — sin 2i • sin 2i' cos (fi - J2 ' ) + 7 * s ^ n i s * n ^ cos 2 (J2 - S7 ') 

lb lb 


(64) 

where primed letters express the orbital elements referred to the equator of the moon 
and the sun, and terms due to the sun and the moon must be added. 

When the longitude of the ascending node of the moon referred to the ecliptic is 
expressed by N, the following relations hold: 
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for the moon, 

sin^ i' = 0. 164 + 0. 068 cos N , 
sin i' cos Q' = 0. 396 + 0. 086 cos N , 
sin i' sin f2' = 0.086 sin N ; 
for the sun, 

sin i' = 0. 399 , 

sin^ i' = 0. 158 , 

n' = 0 . (65) 

If the new variables 


p = sin i cos Q , 

q = sin i sin , (66) 

are adopted instead of i and £1, the equations describing the motion of the orbital plane 
are derived as 


~ |^ = nX 10 5 (5. 06q + 0. 14q + 0. 09 sin N) , 
na q 

^ = -^r - nX 10“ 5 (5.06p - 0. 14p + 0.09 cos N - 0.646) . (67) 

at z dp 
na * 

• _4 

As N = - 1. 47 X 10 n, the solutions of equations (67) are derived by introducing 
two integration constants, -y and t^, as, 
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( 68 ) 


p= 1.06 y cos [-5.06X 10~ 5 n(t - t Q )J + 0. 128 + 0.008 cos N , 
q = y sin [5. 06 X 10 5 n(t - t^)] + 0.008 sin N . 

COS 

If N terms are neglected and y is equal to zero, a stationary solution can be 
obtained: 


n = 0° 

i = T. 3 . (69) 

Other solutions are given approximately by ellipses around the stationary one on 
the (p, q) -plane with a period of 54 years . Therefore, the inclination cannot remain 
very small even if the initial inclination is very small. To keep the inclination below, 
say, one degree for a long interval of time, the initial orbital elements chosen should 
be i = 1° and fi = 270°. 

However, the eccentricity can remain below 0. 001 because a term similar to cos n 
that produces 0. 128 in p comes from the parallactic term with factor e'. 
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